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Abstract. 

Observations of both gamma-ray burst sources and certain classes of active galaxy 
indicate the presence of relativistic shock waves and require the production of high 
energy particles to explain their emission. In this paper we first review the basic 
theory of shock waves in relativistic hydrodynamics and magneto-hydrodynamics, 
emphasising the astrophysically interesting cases. This is followed by an overview of 
the theory of particle acceleration at such shocks. Whereas, for diffusive acceleration 
at non-relativistic shocks, it is the compression ratio which fixes the energetic particle 
spectrum uniquely, acceleration at relativistic shocks is more complicated. In the 
absence of scattering, particles are simply 'compressed' as they pass through the 
shock front. This mechanism - called shock-drift acceleration - enhances the energy 
density in accelerated particles, but does so without changing the spectral index of 
upstream particles. Scattering due to MHD waves leads to multiple encounters between 
the particles and the shock front, producing an energetic particle population which 
depends on the properties of the shock front and the level and nature of particle 
scattering. We describe the method of matching the angular distributions of the 
upstream and downstream distributions at the shock front which leads to predictions of 
the spectral index. Numerical simulation of particle transport provides an alternative 
means of calculating spectral indices, and has recently been extended to cover ultra 
relativistic shocks. We review these calculations and summarise the applications to 
the astrophysics of relativistic jets and fireball models of gamma-ray-bursts. 



PACS numbers: 00.00, 20.00, 42.10 
1. Introduction 

Within the space of thirty minutes, the power emitted in very high energy gamma-rays 
by an active galaxy at a distance of several hundred million light years has been observed 
to increase more than twenty fold ]IJ. This unexpected observation has raised new 
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questions concerning the nature of 'gamma-ray blazars', as this class of source is called. 
In a parallel development, the discovery that at least some of the enigmatic, powerful 
so-called 'gamma-ray burst' sources lie at cosmological distances [0 was perhaps less 
surprising to many astrophysicists, but nonetheless places stringent restrictions on the 
conditions prevalent in the emission regions. These two recent developments have been 
responsible for a renewed focusing of attention on the subject of relativistic shock fronts. 
In each case, the observations imply that the photon density in the source is so large, 
that the gamma-rays would be subject to strong absorption due to photon-photon 
interactions, were not the source moving in our direction with a speed close to that 
of light. The implied Lorentz factors are about 50-100 in the case of the active galaxies, 
and several hundreds in the case of gamma-ray burst sources. 

The subject is not new. Relativistic motion was proposed as a way out of the 
rapid variability problem over thirty years ago ||. The observational motivation then 
was similar to that of today: in order for the variable synchrotron emission in the radio 
band observed from active galaxies to emerge without being subject to self-absorption an 
improbably weak magnetic field was implied. The proposed solution was that the source 
moves at relativistic speed towards the observer, a Lorentz factor of 5-10 being sufficient 
to remove the difficulty. The hypothesis was impressively confirmed by sequences 
of radio observations of the cores of active galaxies, using very long baseline radio 
interferometry, which clearly showed rapid expansion. This apparently 'super luminal' 
motion is now seen in many sources, and indicates the outflow of plasma containing 
relativistic electrons and magnetic fields with bulk Lorentz factors between roughly 5 
and 20 [g. 

Shock fronts are an inevitable consequence when such flows encounter the ambient 
material in the interstellar medium of the host galaxy, or in the intergalactic medium. 
Indeed, they seem the natural sites for the dissipation of the kinetic energy associated 
with the relativistic motion itself. Furthermore, wherever we directly observe shock 
fronts in astrophysics - from the Earth's bow shock to supernova remnants - they 
are associated with particle acceleration. The radiation from both gamma-ray bursts 
and gamma-ray blazars requires the presence of particles (most likely electrons and/or 
positrons) with very high Lorentz factors (exceeding 10 6 in blazars), so that particle 
acceleration must be happening there too. The aim of this review is to summarise 
our current theoretical picture of relativistic shocks in astrophysics and, in particular, 
our ideas about the way in which they might accelerate particles. Section ^| presents an 
overview of hydrodynamic and MHD shocks, Section [3] examines shock drift acceleration 
as well as analytic and numerical methods of treating Fermi acceleration and Section [| 
briefly summarises the applications to relativistic jets and gamma-ray bursts. 
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2. Relativistic Shocks 

Relativistic effects can be important at shock fronts for two distinct reasons. Firstly, the 
post shock temperature can be so high that the thermal motion of individual particles 
approaches the speed of light. Secondly, it is possible that so much energy is released 
into a region containing relatively little matter that the subsequent expansion proceeds 
at a bulk speed approaching that of light. It is certainly possible for the former effect to 
arise when all bulk velocities are quite modest. For example, a hot plasma in which the 
pressure is provided predominantly by photons is always relativistic in this sense, and 
behaves as a gas whose ratio of specific heats is 4/3. If the mass density in such a plasma 
is dominated by atoms or ions, there exists a sound velocity which is determined by the 
pressure of the photons and the inertia of the ions and which can be small compared to 
the speed of light. Downstream of a strong shock front in such a plasma, the bulk speed 
is, of course, less than this sound speed, whereas the upstream fluid speed (provided it is 
still non-relativistic) is equal to seven times the downstream speed. On the other hand, 
it seems unlikely that in a realistic situation a relativistic bulk flow could be thermalised 
by a shock front without producing thermal velocities which are also relativistic. As 
a consequence, the fully relativistic equation of state of a plasma should be used to 
investigate the possible jump conditions. 

Following the solution of the relativistic shock problem in hydrodynamics by 
Taub ||, the question of which types of shock discontinuities are permitted and what 
the jump conditions are in the ideal MHD picture was investigated by de Hoffmann 
& Teller ||. The full relativistic theory was presented by Akhiezer & Polovin and 
an elegant summary of the most important results has been given by Lichnerowicz ||. 
However, the generality of these treatments makes them less accessible to the plasma 
astrophysics community, and several papers have worked out special cases, sometimes 
using particular approximations to the equation of state || [H| 

In view of this, we present in this Section the basic, astrophysically relevant, 
results on both hydrodynamic and MHD shocks in as simple a form as is compatible 
with completeness. In each case a numerical treatment is in general needed to solve 
for the jump conditions, and we provide in the Appendices a detailed description of 
how such an algorithm can be constructed. In the case of hydrodynamic shocks, the 
approximate formulae for strong shocks and for the ultra-relativistic case are presented, 
their derivations being relegated to the Appendix. The special cases which are of most 
interest in MHD are those of weak magnetic field and of ultra relativistic, perpendicular 
shocks. Again, the relevant results are quoted and the derivations presented in the 
Appendix. We also note that ultra-relativistic switch-on shocks are possible only in a 
very limited parameter range. 

The most important property of a shock front for the acceleration problem is its 
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compression ratio, and it has been known for some time [|T^] that increasing the magnetic 
field strength tends to reduce the compression. We present some figures at the end of 
the section which illustrate this effect. 



2.1. Hydrodynamics 

It is useful at this point to recall the theory of relativistic hydrodynamics and the shock 
wave solutions it admits. This will also serve as a means of introducing the notation we 
will use throughout this paper. 

In the absence of external forces and energy sources, the equations of relativistic 
hydrodynamics can be formulated as the vanishing divergence of the stress-energy tensor 
associated with the fluid: 

V M T^ = (1) 

(e.g., Ul3| ). Neglecting dissipative effects, the stress-energy tensor is diagonal in the 
local plasma rest frame, and is given by 

= wu»u v + pg» v (2) 

Here is the four velocity of the fluid (/i = 0, 1,2,3), and g^ v the metric tensor, for 
which we adopt the convention — h++. The scalars w and p are the proper enthalpy 
density and pressure, i.e., those measured in the rest frame of the fluid, in which 
= (1,0,0,0). 

The problem of solving the Rankine-Hugoniot relations to find the jump conditions 
across a shock front requires one to use an equation of state to find the quantity p/p, 
given the quantity e/p. For a dissipation free ideal gas, the Synge equation described in 



Appendix A| is appropriate, and it is necessary to solve Eq. (A9) numerically. Although 
this procedure is straightforward, a popular short-cut is to define a parameter 7 via the 
equation 

p=(7-l)(e-p) (3) 



14, 15], which replaces the equation of state (A9). In the non-relativistic case 7 = 5/3 
and can be identified as the ratio of specific heats of the gas. For a gas whose pressure 
is dominated by a relativistic component, one has 7 = 4/3 (a fully relativistic gas 
has, in addition, e ^> p). For a gas consisting of equal numbers of electrons and 
protons in thermodynamic equilibrium with each other, there is a range of temperatures 
m p ^> T ^> m c in which 7 « 13/9 (here m e is the electron mass and m p the proton 
mass) . 

The jump conditions across a relativistic hydrodynamic shock are given by Eq. (1), 
which expresses the laws of energy and momentum conservation, together with the law 
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of conservation of particle number. For those constituents which are neither created nor 
annihilated in the shock, one has 

V M (n i «") = (4) 

The shock wave is a surface in space-time, = 0, across which there is a 

discontinuity in the fluid variables. The shock normal is given by = d^cf) and, without 
loss of generality, we normalise so that IH^ = 1. Two frames of reference which are 
important from a physical point of view are those in which the upstream and downstream 
plasma is at rest. Defining the Lorentz scalars i>_ and v + to be the shock speed measured 
in the upstream and downstream rest frames respectively, one has 



\ U ±<>H\ /r\ 
V± = I i u ,2 ' ( 5 ) 

where and w+ are the 4-velocities of the plasma immediately upstream and 
immediately downstream of the shock front, respectively. 

Whenever necessary, we will use Cartesian coordinates in which the shock normal 
is along the x-axis. In the upstream and downstream rest frames Z M = (T±v±, F±, 0, 0), 
where r± = (1 — f|) -1 / 2 . As we show below, for a hydrodynamic shock front, the 
downstream plasma flows along the shock normal, as seen from the upstream rest frame, 
so that in this frame w+ = (r re i, — r re ii> re i, 0, 0), where 

V- - V± /r*\ 

Vrel = ~ (6) 

1 — V_V + 

is the relative speed of the upstream gas with respect to the downstream gas, and its 
related Lorentz factor is r re i. 

From the conservation laws, the shock jump conditions are 

M Z„ = (7) 

[TH^ =0- (8) 

where (7) is valid for all conserved constituents. For non-conserved particles additional 
information is required specifying the number density, for example via the equation of 
state. In the shock rest frame, with the shock normal and upstream flow velocity along 
the x-axis it is straightforward to show that the downstream flow velocity is also directed 
along the same axis. In this frame equations (7) and (8) can be written as 

r_p_w_ = T + p + v + (9) 

V 2 _w-vl + v- = r+^+v+ + p+ (10) 
T 2 _w_v_ =T 2 + w + v + (11) 

Given V- and the upstream state, e_/p_, these equations are to be solved for the 
downstream state v+, e + /p_ and the proper compression ratio R = p + /p_. In general, 
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this entails a numerical procedure, which is described in |Appendix A| . However, there 
are several interesting special cases with analytic solutions: 

(i) In the limit of non-relativistic fluid speeds, tij. « 1, the well-known Rankine- 
Hugoniot conditions are recovered. In terms of the upstream Mach number 

/*■ \-l/2 

M_ = «_ | ] (12) 



P- 



one has 



P+ 

v+ 

P± 

P- 



(7 



(7 + 1) 7(7 + l)M 2 _ 

(7-1) 2 
(7 + 1) (7 + 1)M2 

7 + 1 



(13) 



7 - 1 + (2/M£) 

where in our case of a monatomic gas, the relevant value of the adiabatic index is 
7 = 5/3. 

(ii) For a shock in a relativistic gas, in which p = e/3 (both upstream and downstream) 
one finds the simple relation 

v-v+ = - (14) 

However, this corresponds to the perhaps less interesting case of a relativistic shock 
moving into unshocked gas in which the total energy density significantly exceeds 
that corresponding to rest-mass i.e., e_ 3> p_. 

(iii) For a strong shock, at which the upstream pressure can be neglected, the equation 
of state in the upstream gas is unimportant. Furthermore, if particles are conserved 
at the shock, the average Lorentz factor of a particle does not change across the 
shock as seen from the downstream rest frame, i.e., e + = r rel p + (see [Appendix A| , 
(Eq. A20) with 77 = 1). Then, using, for the downstream medium, the equation of 
state for fixed adiabatic index (Eq. 3) one may derive 



™+/Pa 



7 (r rel - 1) + 1 
( w+ / P+ ) 2 (r rel + i) 



7 (2- 7 )(r rel -i) + 2 (15) 

which gives the shock speed and downstream pressure in terms of the relative 
velocity of the upstream and downstream fluids. These relations are generalised in 
[Appendix A| to the case in which particles are created at the shock front. 
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(iv) In the ultra-relativistic limit, T_ — > oo, the upstream pressure may be 

neglected in Eq. (10). If, in addition, the downstream particles are ultra-relativistic, 
in the sense that e + 3> p+, one may combine Eqs. (10), (11) and (6) to find: 

t>+ ->7-l = 1/3 (16) 

r rel - r_^(2-7)/7 = r_/v^ (17) 

These relations are independent of the equations of state upstream and downstream 
and hold whether or not particles are conserved at the shock, provided only that 
the downstream particles are ultra-relativistic. 



2.2. Magnetohydrodynamic Shocks 

In ideal, relativistic MHD, it is assumed that the plasma is dissipation free and that in 
the local rest frame the electric field vanishes. In this case, the electromagnetic field is 
specified by the magnetic field alone, and the system is described by the four 'source- 
free' equations of Maxwell: V ' ll {*F tJ,u ) = where *F^ U is the dual electromagnetic field 
tensor, 






B 1 


B 2 


B 3 \ 


B 1 





E 3 


—E 2 


B 2 


—E 3 





E l 


B 3 


E2 


—Ei 


J 



(18) 



where B = (Bi,B 2 ,B 3 ) and E are the magnetic and electric fields as measured in 
a frame where the plasma's three velocity is f3, augmented by the generalised (ideal) 

— * 

Ohm's Law E = —(3 x B, which constrains the electric field to vanish in the plasma 
rest frame. Defining the four- vector B^ = — u 1/ (*F jUJ/ ) the field tensor can be written as 
^pta> _ B^yU _ u ^b u so that the source-free Maxwell equations become 

V M (BW - u»B v ) = 0. (19) 

In the rest frame of the plasma _B M = (0, B) where B is the magnetic field three vector 
in that frame. In the following, B is taken to denote the magnetic field strength in the 
local plasma rest frame which satisfies B^B^ = B 2 . The components of B^ in a frame 
where the plasma is moving with four velocity can be derived from the appropriate 
Lorentz transformation in terms of B. The energy momentum tensor of the system 
consisting of electromagnetic fields and fluid is 

/ B 2 \ ( B 2 \ B^B U 

r-=(. + _)„.v+( P + _) 9 --_. (20) 

and the equations of relativistic MHD consist of the vanishing divergence of this tensor, 
together with (19) and, of course, an equation of state for the fluid. An important 
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difference between relativistic and non-relativistic MHD is that in relativistic MHD it is 
not possible simplify (19) by ignoring the displacement current. This is because in the 
relativistic case, fluctuations in the space charge cannot be neutralised on a time scale 
which is much faster than the other dynamical time-scales. 

There are three wave modes in MHD, the fast and slow magnetosonic modes, and the 
Alfven mode. The phase velocities of these modes in a given direction will be denoted 
by ffast, v s \ ow and v&. Quite generally one has 



Vfast >v A > v s{ow 



(21) 



Defining $ as the angle between the magnetic field and the direction of propagation, 
measured in the plasma rest frame, the Alfven wave has the speed 



va 



E 2 /(4tt) 



cos $ 



\\ w + B 2 /(4ir) 

and the fast and slow mode speeds are the roots of the equation 



w(v s 2 — l)f 4 — 



w + 



B 2 



4irv 2 



B 2 

l-v 2 ) + — cos 2 $(l -w 2 ) = 

477 



(22) 



(23) 



|i0| . These wave speeds are important for the properties of shock fronts [0, §| , since, as 
in the non-relativistic case Jl6|, |17|? a physically realisable shock front cannot provide 
a transition across the Alfvenic point in ideal MHD. In other words, the two kinds of 
transition available are the slow-mode shock, with 



Va- > V- > U s iow- 
> ^slow+ > V + 

and the fast-mode shock which has 

V- > Ufest- > ^A- 

v ias t+ > v+ > v A + 



(24) 



(25) 



where v f as t,A,siow± refer to the propagation speeds along the normal to the shock front in 
the up and downstream regions. The fast-mode shock is associated with an increase in 
the magnetic field strength across the shock, whereas the slow-mode shock results in a 
decrease of the field. 

For the purposes of numerical simulation of a relativistic MHD flow, it is interesting 
to put the governing equations into divergence form, avoiding the additional algebraic 
constraint implied by the relation u^B^ = 0. Such a formulation has been presented 
T8| and implemented [l£j by van Putten. However, at least for the non-relativistic case, 



other methods are available (e.g., Falle et al. |[20|| ) and there is also no advantage to be 
gained by such a reformulation in the investigation of the jump conditions. These are 
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determined by the equations of conservation of particles and of energy /momentum (7) 
and (8) [as in the hydrodynamic case, but with the modified stress-energy tensor of 
Eq. (20)], supplemented by equation (19) which implies: 

[B»u u - u»B v ] 1^ = 0. (26) 

It proves convenient to define 5± = B±/ v / 47rp_ so that \b\± = B±/4np_ is twice the 
ratio of the magnetic field energy density in the local plasma rest frame to the upstream 
rest mass energy density. In the upstream and downstream rest frames we define $± to 
be the angle between the magnetic field and the shock front's direction of propagation, 
and take coordinates such that the field lies in the x-z plane. Therefore we can write 
for the two vectors b + and b_ 

= (47rp_)- 1/2 (0, B± cos $±, 0, B± sin $ ± ) (27) 

in the upstream or downstream frame as appropriate. Given t>_, $_ and \b\ 2 _, the jump 
conditions can be solved numerically to find v + , $+, \b\^_ and the compression ratio 
p+/p_ using the method of Majorana & Anile fllQfl , which is described in |Appendix E| , 



and involves a single one- dimensional root-finding operation. 

In general, the plane which contains the shock normal and the upstream magnetic 
field, as seen in the upstream rest frame, also contains the plasma velocity and 
the magnetic field in the downstream half-space - i.e., these shocks do not permit 
the generation of a non-coplanar magnetic field. However, in addition to the fast- 
mode and slow-mode shocks, there exist solutions of the jump conditions (known as 
rotational discontinuities), which propagate with the Alfven speed and do not involve 
a compression of the plasma, or a change in the field strength, but merely a change 
in the direction of the field and the plasma speed. In the standard picture of Fermi 
acceleration, they are not expected to accelerate particles. Slow mode shocks are also 
thought to be ineffective accelerators. Not only does the magnitude of the magnetic 
field decrease across the shock, (rendering shock drift acceleration ineffective) but also, 
if particle scattering is produced by Alfven waves, then the speed of the scattering centres 
in the fluid frame exceeds the speed of the slow mode shock. Crossing and recrossing 
would then lead to a net deceleration if the waves predominantly stream away from 
the shock front, which would turn diffusive acceleration into diffusive deceleration |2l[ . 
Therefore, only fast-mode shocks are considered in the following. 

As in the hydrodynamical case, there are some interesting limiting cases for MHD 
shocks: 

(i) In the case of a weak, dynamically unimportant magnetic field, such as considered 
by Kirk & Heavens p2 |, Begelman & Kirk [23] several useful analytic results can 



be obtained. Of course, the jump conditions for the plasma speed and temperature 
are identical to the hydrodynamical case. In the shock rest frame with both v_ 
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and v + directed along the x-axis, the x-component of magnetic field is conserved 
across the shock. Since the component of magnetic field in the direction of the flow 
velocity is unchanged by a Lorentz boost to the local fluid frame this gives a simple 
expression connecting the angles and $+, 

B_ cos $_ = B + cos $+ (28) 

The (proper) compression in the magnetic field can be related to the upstream 
parameters and the compression ratio r = v_/v + = RT + /T_ of the fluid: 



(ii) If the upstream plasma flows along the magnetic field and the shock normal, it is 
possible to have a purely hydrodynamic transition in which the downstream flow 
is also along the magnetic field. However, for some parameter ranges, it is also 
possible to generate a transverse component of the field at the shock front in what 
is called a 'switch-on' shock. This happens if the hydrodynamic solution violates 
the constraints given by (24) and (25) (referred to as evolutionarity conditionsf). 
To find the relevant range of parameters in the ultra-relativistic case, we make the 
approximations 

V- — > 1; v+ — > 1/3; r+ -> 3%/2/4 

T rel -> T_/V2; w+/p + -> 2v / 2T'_/3 (30) 
[using Eq. (15]. For a parallel shock, where 6_ = b + , the Alfvenic Mach number is 

m " = t(^ + " 2 ) 1/2 (31) 

A switch-on shock occurs if the upstream flow is super-Alfvenic and the 
hydrodynamic jump conditions give a sub- Alfvenic downstream speed, which 
means, in the ultrarelativistic approximation, 

Vl + b 2 <T_<V3b (32) 

This parameter range is unusual, since it requires a plasma in which the magnetic 
field energy density greatly exceeds the rest mass energy density, as well as a shock 
front whose normal is closely aligned with the upstream magnetic field. 

f Lichnerowicz has pointed out that also the switch-on (and the switch-off) shock is not 
evolutionary. The physical reason for this is that the plasma cannot break the azimuthal symmetry 
around the direction of the magnetic field in order to create the transverse component. However, for 
an arbitrarily small angle between the shock normal and magnetic field, this symmetry is broken, and 
a shock appears which is both evolutionary and arbitrarily 'close' to the switch-on shock. On the other 
hand, the purely hydrodynamic solution does not have an evolutionary shock in its neighbourhood, 
when conditions (24) and (25) are violated JlCfl . 
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(iii) The ultra-relativistic approximation is more commonly encountered in connection 
with perpendicular shocks, such as are thought to occur in the winds driven by 
radio pulsars such as the Crab [25, [12], 2j| ^7]. Here it is convenient to use the 
(Lorentz invariant) parameter 

B 2 



a 



Attw 

= p^b 2 /w (33) 

which, in the case of a magnetic field orientated perpendicular to the flow, is just 
the ratio of the energy flux carried by the electromagnetic field to that carried by 
the particles. Kennel & Coroniti |2£| have given an approximate solution for strong, 
ultrarelativistic perpendicular shocks: 



8a± + 26a_ + 17 + J64a± (a 



1) 



20ct_((t_ 



D + l 



p+ 



p+ 

p- 

1 + 



r 2 



i + (j- i 



16(<r- 



i) 



(34) 



which goes over into the hydrodynamic solution (15) and (29) for er_ — > 0. 
As pointed out by Kundt & Krotscheck |T2|], the compression which can be obtained 



in a relativistic shock decreases as the magnetic field becomes more and more important 
dynamically. To illustrate this point we present in figure [I] the jump conditions solved 
for the case of a magnetic field at an angle of 45° to the shock normal in the upstream 
frame. The shock is taken to be strong, in the sense that the pressure of the plasma 
upstream is negligible. However, we allow for an appreciable magnetic pressure in the 
upstream plasma, as measured by the Alfvenic Mach number, defined as Ma = t>_/t>A, 
with the (relativistic) Alfven speed given by equation (22) . The Synge equation of state 
is used in these computations, and at Alfven Mach numbers above 100, the result is 
close to the hydrodynamical case. In addition, figure [l] shows the hydrodynamic result 
for a strong shock using the fixed adiabatic index approximation of equation (15), and 
taking the ultrarelativistic value 7 = 4/3 (which, of course deviates significantly from 
the Synge equation of state in the non-relativistic case). 

These computations were performed using an algorithm based on the description 
in |Appendix B| . They extend the results presented by Ballard & Heavens [^] by using 
the Synge equation of state and showing explicitly that the division between super- and 
subluminal shocks used in that paper (and discussed in Section |3|) is unimportant for the 
jump conditions. Similar results, but in terms of the parameter <r_ = v 2 _/(Mj L _ — v 2 _), 
are given for the ultra-relativistic perpendicular shock by Kennel & Coroniti 
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Figure 1. The compression ratio of a strong oblique fast-mode shock front ($_ = 45°) 
propagating into a magnetised plasma. Various Alfven Mach numbers Ma are shown, 
together with the hydrodynamic approximation of equation (15) (dotted line). 



3. Particle acceleration 

Restricting ourselves to the picture of a shock front as a discontinuity in an MHD 
flow, the crucial question to be answered before particle acceleration can be discussed 
is the precise nature of the interaction between accelerated particles, the plasma and 
its embedded magnetic field. The standard starting point is orbit theory. Provided the 
embedded magnetic field varies slowly in time and space compared to the gyro-frequency 
and gyro-radius of the particle, the orbit can be described by the location of the gyro- 
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centre (three spatial coordinates) the energy, the first adiabatic invariant or magnetic 
moment and the gyro-phase. The latter, which is the angle variable corresponding to the 
action variable 'magnetic moment', varies on the fastest time-scale. Correspondingly, 
it is assumed that the distribution function quickly adjusts itself to be independent of 
this coordinate. From this 'unperturbed' orbit, the effect of small fluctuations in the 
electromagnetic field can be taken into account using quasi-linear theory. The result 
is a transport equation which is a second order partial differential equation describing 
diffusion in energy, magnetic moment (= p(l — [i 2 )/B, where p is the momentum and [i 
the cosine of the 'pitch angle' between the velocity vector and the magnetic field) and 
a coordinate perpendicular to the magnetic field. Explicit expressions for the diffusion 
coefficients, which are applicable in the rest frame of the plasma, have been computed 
by several authors [[29|, 0, [H], |3^| . The most important wave mode responsible for the 
fluctuations is thought to be the Alfven wave propagating parallel to the magnetic field, 
since it is only weakly damped |33|, |34j. In this case the fastest diffusion time-scale is 
that associated with pitch angle scattering. Terms involving diffusion in energy and in 
position of the gyro-centre are smaller. In the next step, it is assumed that pitch angle 
diffusion is rapid enough to keep the particle distribution function close to isotropy, 
which enables the transport equation to be reduced to one containing only diffusion in 
the spatial coordinate along the field line, together with smaller terms describing cross- 



field and energy diffusion [J35], [36|]. Neglecting diffusion in energy, one arrives at the 
well-known cosmic ray transport equation [37|, [38]]. Applying it to a flow pattern which 
contains a non-relativistic shock front leads to the theory of diffusive shock acceleration 
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This acceleration mechanism has been intensively studied (for reviews see Drury [f£| 
and Kirk et al. [EH]). At a non-relativistic shock wave for which t>_ and v+ are each very 
much less than the particle speed v p , the particle distribution at the shock is, to lowest 
order in v±/v p , isotropic when viewed in either the upstream or downstream frames. A 
particle which crosses from upstream to downstream and back again then undergoes a 
small fractional increase in its momentum which is proportional to v m \/ v p and is given by 
Ap/p = 4(t>_ — v + )/3v p . The factor 4/3 is a consequence of isotropy. The only process 
which can stop this slow but inexorable acceleration is escape of the particle from the 
vicinity of the shock. With n s the number density of particles at the shock the flux 
of (isotropic) particles of speed v p crossing the shock is n s v p /A. Diffusion and isotropy 
guarantee that the number density of particles is constant in the downstream region so 
that the flux of escaping particles is then n s v + . The probability of escape per shock 
crossing is then the ratio of the flux of escaping particles to those crossing the shock, 
Pesc = 4f + /fp. The ratio of P esc to Ap/p determines the slope of the integral particle 
spectrum (i.e. the number of particles above a momentum p). The phase space density 
of particles then has a power law spectrum f(p) oc p~ q where q = —P esc /(Ap/p) — 3. 
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This gives a spectral index, which depends on the compression ratio of the shock alone, 
of q = 3r/(r — 1) where r = v^/v + . Note that this result does not depend on the level 
of scattering as long as that scattering is diffusive leading to an isotropic distribution 
when seen in each fluid frame. 

Relativistic flow brings substantial changes to this picture, because the assumptions 
made in deriving the spatial diffusion equation are no longer valid. Thus, when v m \ is of 
order the velocity of light, the quantity v Te \/v p can no longer be small. The Lorentz boost 
linking the two frames significantly deforms the angular distribution of the particles, so 
that it cannot be nearly isotropic in both the upstream and downstream rest frames. 
Instead of assuming spatial diffusion, the angular dependence of the distribution function 
must be computed explicitly. 

Nevertheless, one can gain some understanding by exploiting a key property of the 
quasi-linear approach - namely that diffusion across the field lines is ineffective - despite 
the fact that a firm theoretical basis is lacking in the relativistic case. This property 
implies that, to a first approximation, a particle's guiding centre remains fixed to a 
given field line except at points where its orbit cuts the shock front. The kinematic 
implications are severe. The pattern of successive crossings and re-crossings of the 
shock front as a result of diffusion along the field lines, becomes progressively more 
difficult to maintain as the obliquity of the shock increases, since a diffusing particle 
has to chase not the shock front, but the point of intersection of the shock front with 
the field line to which it is tied. Re-crossing becomes impossible when the intersection 
point moves with a speed greater than that of light, i.e., t> int = t>_/ cos$_ > 1, and 
this divides shock fronts into two categories: subluminal (vi nt < 1) and superluminal 
(f int > 1). The jump conditions themselves (see Section do not exhibit any special 
features at the change-over between sub- and superluminal shocks, however, there are 
interesting differences when one considers the properties of the two types of shock under 
Lorentz transformations. Thus, for subluminal shocks it is always possible to find a 
transformation to the de Hoffmann/Teller frame in which the electric field vanishes 
everywhere. This is a convenient frame for computational purposes. In the case of 
superluminal shocks, no such simple frame can be found. The best that can be done 
is to transform to a frame in which the magnetic field is perpendicular to the shock 
normal, the shock is stationary, and the incoming flow is in the same plane as the shock 
normal and the magnetic field |43|, ^3|. Generally, if one assumes the direction of the 
magnetic field in the upstream plasma has no causal connection with the shock front, 
so that it is randomly oriented with respect to the shock normal, subluminal shocks 
will greatly outnumber superluminal ones in subrelativistic flows. On the other hand, 
in relativistic flows subluminal shocks will be rare, since they require orientation of the 
magnetic field to within an angle of 1/T_ of the shock normal. 

At superluminal shocks, perturbations to the particle orbit are unimportant for 
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acceleration when cross-field transport is neglected. However, even in the absence of 
fluctuations which scatter energetic particles, a shock front in a perfect fluid is capable 
of increasing the average energy of a distribution of particles incident upon it. Seen 
from the upstream rest frame, there exists an electric field in the downstream half-space 
simply because it is in motion (the — v x B field). A charged test particle can tap into 
this energy source if its gyrations through the regions of different magnetic field strength 
upstream and downstream cause it to drift along the shock front in the appropriate 
direction. This is usually called 'shock-drift' acceleration. It can be analysed by tracing 
a large number of orbits numerically. The relativistic version of this mechanism has 



been investigated in |23[ and is described in Section 3.1 



In Section [3]2] we turn to the situation in which scattering does play an important 
role in the acceleration mechanism - i.e., the relativistic generalisation of the 
diffusive shock acceleration mechanism - and summarise the analytic and semi-analytic 
approaches to the problem. These involve the solution of the angular dependent 
transport equation assuming an explicit form for the pitch angle scattering coefficient. 
As such, they neglect cross field transport and apply principally to subluminal shocks. 
However, these methods may also be applicable to the case in which no a priori 
assumption is made concerning cross-field diffusion. If the direction of the embedded 
magnetic field is unimportant because a particle performs only a fraction of gyration 
between shock crossings (as in the case of an ultra- relativistic shock), or because the 
fluctuations in the field strength reach or exceed the amplitude of the embedded field, 
then transport by successive small angle deflections can be treated simply by replacing 
the pitch angle by the angle between the shock normal and the particle velocity |45[ . 

Despite the advantages of a semi-analytic approach, Monte-Carlo simulation 
provides the most popular technique for investigating acceleration at relativistic shocks. 
The formulation of the problem is conceptually simple, and the method permits one to 
look at a variety of effects (such as synchrotron losses and time dependence) which have 
so far defied an analytic approach. This work is discussed in Section |3.3j . 

The ideal MHD approach, coupled with a simple transport prescription for energetic 
particles allows us to make some predictions - such as spectral indices - which are, 
at least in principle, accessible to observation. Currently, the only work of which 
we are aware which goes beyond this picture concerns the structure of relativistic 
perpendicular shocks, with particular application to the Crab Nebula j|6], [47], pE8| . The 
techniques involved - one and two dimensional particle-in-cell or hybrid simulations - 
differ markedly from those employed in the perfect fluid case, and are not discussed in 
this review. 
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3.1. Shock-drift acceleration 

The dramatic increase in surface brightness which can be produced by a relativistic 
shock front merely as a result of the 'compression' of the electrons was pointed out 
by P. Scheuer |[49|| . This raises the possibility that observations of hot-spots in extra- 
galactic radio sources could be understood even if the Fermi process with its emphasis 
on crossing and recrossing of the shock front were unimportant. 

To quantify this statement, consider a gas of relativistic electrons with an isotropic 
distribution function in the local fluid frame. They achieve isotropy by experiencing 
elastic scattering by slowly moving, low-frequency MHD waves, which may be self- 
excited. The distribution in Lorentz factor 7 is unaffected by the scatterings. 

The phase space density of particles contained in the interval dp around p is 
conventionally assumed to be a power law: 

f(p) = Cp~ s (35) 

between a lower and an upper cut off: p min < p < p nmx . The corresponding number 
density within this range of p is n(p) = 4irp 2 Cp~ s . Let us assume that this distribution 
with C = C- accurately describes particles in a fluid element which flows into an MHD 
shock. 

If the scattering events are so rapid that the length scale over which they isotropise 
the electrons is much shorter than the length scale characterising the thickness of 
the shock, then the relativistic electrons react adiabatically. The Fermi process is 
unimportant in this case because scattering anchors the particles in the local fluid 
element and prevents them from repeatedly sampling the velocity difference between the 
upstream and downstream sides of the shock in between scatterings. The characteristics 
of the transport equation are then 

pp^ 1 ^ 3 = constant (36) 



where p is the proper fluid density (e.g., pip, so that downstream of the shock front 
one has 

U(p) = f-[p(p+/p-) 1/3 ] 

y\ P- S (37) 

Provided p is far from the cut-off momenta, this can be written in terms of a (proper) 
compression ratio for the electron distribution R e : 

f-(p) 

= R s/3 (39) 
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where R is the proper compression ratio of the MHD fluid. 

Even if they are not scattered whilst traversing the shock, charged particles may 
still be said to behave 'adiabatically' in the sense that the first adiabatic invariant of 
motion (the magnetic moment) is conserved. For this approximation to be valid, it 
suffices that the length scale on which the magnetic field and fluid speed change (i.e., 
the shock thickness) is long compared to the gyro-radius of the electron concerned. 

Perhaps surprisingly, the magnetic moment is still conserved to a good 
approximation even if the magnetic field changes discontinuously. This was pointed out 
by Parker |50| and Schatzman |5l| for perpendicular shocks and investigated numerically 
by Decker |52 |. From the work of Whipple et al |53] it is known that the magnetic 
moment is conserved during encounter with a discontinuous change in the magnetic field 
at a perpendicular shock front provided the velocity of the front is small compared to that 
of the electron. In the case of synchrotron emitting electrons, this means d_<1. Away 
from the shock front (i.e., further than a gyro-radius from it) the particle distribution 
may relax slowly to isotropy as a result of pitch angle scattering. However, crossing and 
re-crossing is ruled out if the shock is superluminal, in which case this process operates 
undisturbed by the Fermi mechanism. The effective compression ratio for the electron 
distribution defined in equation (38) is straightforward to calculate using Liouville's 
theorem applied to the drift-kinetic equation (e.g., 



R. = Ii$-W Q -WbJI^) (40) 



2 a V2' 2 

7T ( ._ 1)/2 r(( 3 -i)/2) 

2 B r(s/2) 
for R B > 3 (41) 

where Rb is the proper compression ratio of the magnetic field, to = (Rb — 1)/Rb and 
B Q (a, b) is an incomplete Beta function (f54|, page 944). The special case s = 4 can be 
expressed in terms of elementary functions, and is given by van der Laan [55fl. If the 
magnetic field is dynamically unimportant, the magnetic compression Rb is given by 
(29). 

At a relativistic shock, the magnetic moment of a charged particle is in general 
not even approximately conserved, although for subluminal shocks there still exists a 
small region of phase space for particles moving with the shock front within which the 
approximation holds | 56| . Notwithstanding this, significant progress can only be made 
by numerically following individual particle trajectories as they are carried over the 
shock front. Since those sections of the path between intersections are simply helices, 
the algorithm reduces to locating the intersection points and performing the relevant 
Lorentz transformations. Phase is now an important coordinate, so that results must 
be integrated over all initial phases, as well as initial pitch angles. The results presented 
in p5) show that the electron compression ratio R e is much larger at relativistic shocks 
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than a naive extrapolation of the non-relativistic results (39) and (41) would indicate. 
In fact, for r rc i 1, one has R c oc T^ 2 , compared to the asymptotic behaviours 
R c oc r*gf and R c oc predicted by equation (39) and equation (41), respectively. 

In connection with this result, it is interesting to note that at a strong relativistic shock 
front, the average Lorentz factor of thermal particles downstream is the same as their 
pre-shock value, as seen from the downstream rest frame. Essentially, the energisation 
in shock drift acceleration at an ultra relativistic shock front is only slightly greater than 
is achieved by isotropising the upstream particles using the downstream magnetic field 
and the static fluctuations associated with it. 

Although each of the three mechanisms described in this section results in a strong 
enhancement of the energy density in accelerated particles, the boost in the energy 
of each individual particle is modest. The mechanisms work only when a population 
of particles with a non-thermal distribution flows in towards a shock front, and the 
characteristics, such as spectral index, of the incoming distribution are not changed. 
Thus, there is no characteristic spectral signature, nor is it possible by these mechanisms 
to accelerate a particle over several decades in energy. Nevertheless, these mechanisms 
would be the simplest ones to incorporate into a numerical hydrodynamic or MHD 
simulation. 

3.2. First order Fermi acceleration - analytic and semi- analytic methods 

In contrast to the shock-drift mechanism, the repeated crossings and re-crossings 
characteristic of first order Fermi acceleration enable large boosts in particle energy 
and, perhaps more importantly, produce a particle distribution with a characteristic 
spectral index. Early on, it was thought that the increased escape probability in the 
relativistic case would mean that most of the accelerated particles would be swept away 
from the shock never to return, which would severely reduce the acceleration 'efficiency' 
(i.e., cause the spectrum to fall off rapidly to higher energies). However, this does not 
happen for two reasons. Firstly, the escape probability is related to the speed of the 
shock in the downstream rest frame. This quantity does not increase indefinitely as the 
shock becomes faster, but remains lower than the relevant wave speed. In the case of 
a hydrodynamic shock, the wave speed is that of sound, which cannot exceed c/y/3. 
Secondly, the acceleration per shock crossing also increases as the shock becomes faster. 
In fact, there is no upper limit to this quantity, since a particle which encounters a 
shock head-on (i.e., whilst moving along the shock normal) and returns in the opposite 
direction suffers a fractional increase in energy which is of the order of r 2 el . However, 
the average fractional energy gain is a sensitive function of the angular distribution of 
the particles, and, as we shall see, this adjusts itself such that the competition between 
acceleration and escape results in a power-law spectrum which does not differ radically 



19 



from the non-relativistic value. 

To find the spectrum of particles accelerated by the first order Fermi process at 
a shock front it is necessary to solve the transport equation in the upstream and 
downstream media and match the solutions at the shock front. In contrast with the 
extensive analytic developments in the non-relativistic theory, only a few special cases 
have been solved for relativistic flows. The simplest of these, which is outlined below is 
the determination of the characteristic power-law index of test particles accelerated at 
a plane parallel shock front in which the fluid flows along the shock normal. This was 
first attacked by Peacock ]57|, who, however, did not solve the full transport problem, 
but assumed a particular angular distribution of particles in the downstream plasma. 
It was subsequently solved using an eigenvalue expansion technique J58|. This method 
has been generalised to treat acceleration at oblique shocks, employing the additional 
assumptions of conservation of magnetic moment upon crossing the shock, and neglect 



of cross-field transport processes [22|. Much harder spectral indices were found in these 



cases, largely as a result of the accumulation of particles swept up into a precursor by 



repeated reflections from the shock front |59| , |60| , |56| . These are ameliorated somewhat 
when account is taken of the weakening of the compression ratio of an oblique shock 
front by magnetic pressure [28|, as shown in figure [|. Steeper spectra may also result if 
cross-field transport becomes important. 

'Modified' (parallel) shocks - i.e., those of finite thickness and a prescribed velocity 
profile formed the subject of another generalisation of the method [^T, |62|]. The 
technique described in the latter paper is relevant for the problem of injection at non- 
relativistic shocks, where the anisotropy of the particle distribution plays a key role 
p3fl . Adding a term describing large-angle scattering to the transport equation allows 
further insight into the acceleration mechanism and its sensitivity to the scattering 
operator using essentially the same technique [p4| . Non-linear effects, which play a 



prominent role in non-relativistic theory, have proved difficult to treat analytically - 
the only advance being the relativistic generalisation |H| of the classification of possible 
stationary solutions given in the diffusive case by Drury and Volk [j36|j . 

Returning to the problem of test particles at a parallel shock, in the presence of 
pitch angle diffusion described by the coefficient D w , but neglecting as usual diffusion 
in energy, the equation to be solved is 

r ± d + ^)f + r±( „ ± + ,)f ^Azvg, (42) 

where we have assumed the Lorentz factor of the particles is much larger than that of 
the flow and have accordingly replaced their velocity with the speed of light (= 1). The 
relatively simple form of (42) is a consequence of a mixed coordinate system |J7| in 
which the cosine of the pitch angle \x is measured in the local rest frame of the plasma, 
but the space-time coordinates x and t refer to the rest frame of the shock front. Seeking 
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-p^—qrM = A^(v + n)QTM ■ (44) 



a stationary solution in the shock frame, we separate the variables /i and x to find an 
expression for the general solution in terms of an eigenfunction expansion: 

oo 

/ = E 9i(p)Q v MeM^/n , (43) 

i=— oo 

where the <7i(p) are arbitrary functions of momentum. This solution is valid in both the 
upstream (x < 0, v = t>_) or downstream (x > 0, v = v + ) regions. The eigenvalues A" 
and eigenfunctions Qi(fi) appropriate to each half space are solutions of the equation 

dfi w d/jb 

together with the boundary conditions that the eigenfunctions be regular at the singular 
points \i = ±1 |58| . 

The eigenvalue problem equation (44) resembles the standard Sturm-Liouville 
problem familiar in many branches of physics and shares with it several convenient 
properties ([0, page 227), such as a discrete spectrum with eigenvalues ordered 
according to the number of nodes displayed by the corresponding eigenfunction. 
However, the special feature of equation (44) is that precisely at that value of // at 
which particles are stationary with respect to the shock front, the 'weighting function' 
v± + fi changes sign. As a result, there are two families of eigenvalues and eigenfunctions 
for each value of v: one of them has positive, the other negative eigenvalues A". It is 
convenient to use the label i to indicate this, choosing positive integers for the positive 
eigenvalues and negative for the negative ones. The number of nodes possessed by the 
eigenfunctions is then |i| — 1. In addition to the two families of eigenvalues, there exists 
a special isotropic eigenfunction with zero eigenvalue which is labelled by i = 0: Aq = 0, 
Qoip) =constant. This just reflects the property that pitch angle diffusion tries to drive 
any distribution towards an isotropic one, which is the only stationary solution of the 
transport equation in the absence of boundary conditions at finite x. The eigenfunctions 
are orthogonal and may be normalised such that 

" +1 d/iQ^ + yU )QJ = % . (45) 

where 77^ = for i 7^ j, rja = 1 for i < and rju = —1 for i > 0. 

In general, it is necessary to use a numerical technique to solve equation (44). Both 
a Galerkin method 0, 69] and a shooting method |7(J have been employed. In the 
limiting cases of non-relativistic flow t)<l and ultrarelativistic flow v — > 1, approximate 
analytic expressions are available, which are described in [Appendix C 



In the presence of a shock front at x — 0, the general solution equation (43) indicates 
that, in the upstream half-space {x < 0), all terms containing eigenfunctions with i < 
diverge at large distance from the shock (x — > —00). Thus, for a physically admissible 
solution, the coefficients g^(p) must vanish for i < 0. Furthermore, if we are interested 
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in particles accelerated by the shock front, we can demand that the density vanish 
far upstream, which implies that the coefficient with % = must also vanish. Similar 
arguments can be applied to the downstream distribution, which, however, has non- 
vanishing density at x — > oo. 

To find an approximation to the characteristic power-law index, one assumes 
/ oc p~ s , so that in equation (43) giip) = p~ s , and uses as an Ansatz an approximation 
to the downstream distribution which fulfills the boundary conditions at x — > oo. 

o 

f(p,p,x) = p~ s E OiQ+(/i)exp(A+x/r + ), (46) 
%=-j 

where the a« are constants, and we have truncated the expansion to include only a finite 
number (J + 1) eigenfunctions. We have also used the notation that the momentum 
and cosine of the pitch-angle measured in the downstream rest frame are (p, //). They 
are related to those measured in the upstream frame (p, /i) by a Lorentz transformation 

p = r rel p(l + v vel fx) (47) 
V = T- , (48) 

1 + V le ifJ, 

The distribution function expressed as a function of downstream quantities is denoted 
by /. On crossing the shock, particles are assumed not to undergo a sudden change in 
momentum (or direction of motion) so that Liouville's theorem demands continuity of 
the distribution. This condition reads 

f(p,fji,x = 0) = f(p,ji,x = 0) (49) 

Accordingly, the upstream distribution at the shock front, x — 0, may be computed 
from equation (46) by substituting for p and p, using equations (47) and (48). One must 
now try to fit the upstream boundary conditions as well as possible using the J + 1 
constants and the power-law index s. Although / of equation (46) cannot be required 
to be orthogonal to all the upstream eigenfunctions which cause divergent behaviour at 
x — > — oo, it can at least be made orthogonal to those J of them with the longest range 
(i.e., smallest |A|). Projecting onto these gives: 

o r+l 

r rei zl / dfiQ i {fi){v^ + fi){l+v rclf i) s Q+(p)aj 

= (50) 

for % = — J...0. The condition that the set of J + 1 homogeneous equations (50) 
possess a non-trivial solution is sufficient to determine the unknown power-law index 
s. Once this is determined, the coefficients dj and hence the angular dependence of the 
distribution follow. 

An example of the angular distribution at a relativistic shock, as seen from the 
rest frame of the downstream plasma, is shown in Fig. 0. This figure was computed 
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Figure 2. The pitch-angle distribution of accelerated particles at a parallel relativistic 
shock front with V- = 0.9, and v + = 0.37, as a function of the cosine /j, of the pitch 
angle measured in the rest frame of the downstream plasma. Isotropic pitch-angle 
diffusion is assumed and the normalisation of / is arbitrary. The depletion of particles 
with /i « 1, (those which move almost along the shock normal into the downstream 
plasma) arises because particles which move into the upstream plasma are overtaken 



again by the shock before undergoing substantial deflection 45 
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using an isotropic pitch- angle diffusion coefficient fl w oc 1 — fi . As well as confirming 
that the the distribution is strongly pitch-angle dependent, Fig. [| shows that very few 
particles travel in the direction /i = 1, i.e., along the shock normal into the downstream 
region. The reason is that a particle which crosses into the upstream plasma undergoes 
relatively little deflection before being caught again by the relativistically moving shock. 

Results for the power-law index s obtained by Kirk Q are summarised in Fig. (|3]), 
where the characteristic spectral index is plotted as a function of t>_r_ for parallel shocks 
having various equations of state. The plasma composition in this case is approximately 
primordial, with 25% helium by mass. Heavens and Drury |70| have given approximate 



formulae for s which they found by fitting their results in the range < v _ < 0.98. 
For isotropic pitch angle diffusion and a gas consisting of electrons and protons in 
thermodynamic equilibrium, they find for a strong shock 

s « 3.99 - 3.16w_ + 10.86t» 2 - Ihvl + 7.46^1 (51) 

This fit is also shown in Fig. The agreement between the two treatments is clearly 
very good; the minor differences probably being due to slight variations in the jump 
conditions brought about by the different plasma compositions. 

The main features of the plot can be understood as simply being due to the 
different compression ratios which arise when the shock fronts become relativistic. 
Nevertheless, significant deviations from a naive extrapolation of the non-relativistic 
formula s = 3r/(r — 1) are found already at quite low speeds, as is shown in Fig. (H). 

In contrast to the diffusive case, the value of the spectral index s for relativistic 
shocks depends on the functional form of the pitch-angle diffusion coefficient D^. This 
point has been investigated by Heavens and Drury 0] and by Kirk ||69|| . For example, 



if pitch-angle scattering through the point // = is severely restricted [jflU , the spectrum 



is steepened. Figure [| illustrates this for pitch-angle scattering given by 

D m oc (1 — /i 2 )/i 9 for > e 

= constant for < e (52) 

with e = 1/30 and the index q, which corresponds to the power-law spectrum of the 
turbulent wave-energy in the quasi-linear theory equal, taken in this example to be 2. 
Heavens and Drury, on the other hand adopt the prescription 

D^ = (l-^X^ + O.Ol) 1 / 3 (53) 

which roughly corresponds to the quasi-linear result in the presence of Kolmogorov 
turbulence. 

To date, the method has been applied to shocks moving with a maximum Lorentz 
factor T_ = 5. Although the results in Fig. [3] seem to indicate a convergence to a value 
around 4.2 - a result also found and commented upon by Heavens and Drury [[70 
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Figure 3. The spectral index as determined using the Qj method for a parallel 
relativistic shock in plasma with various equations of state. Plotted on the abscissa 
is the x component of the four velocity of the upstream plasma v_r_. The dotted 
line corresponds to a relativistic gas with jump conditions given by equation (14). 
The dashed line depicts a gas consisting of fully ionised hydrogen and helium (25% by 
mass) with the electrons assumed cold (T e = 0); the solid line the same gas but with 
cold ions (Tj = 0, i.e., dominated by electron pressure). Two sets of results are shown 
for complete thermodynamic equilibrium (T c = T\), those of Kirk ^] (dashed-dotted 
line) for the same composition as the other cases and those of Heavens & Drury |7(| 
for a proton/electron gas (long dashed line), as given by equation (51). 
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Figure 4. The spectral index as determined using the Qj method (solid line) for a 
parallel relativistic shock in plasma (25% helium) in full thermodynamic equilibrium. 
The dashed line corresponds to the non-relativistic formula for the spectral index: 
s = 3r/(r — 1), where r = v-/v+ is the compression ratio 



there is no analytic guarantee that the asymptotic limit either exists or is approached 
smoothly. Notwithstanding this, recent numerical results ( [72], f|5| - see section |Q| ) also 
find convergence to s = 4.2 for very large Lorentz factors. 



3.3. Numerical treatments of acceleration 
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Figure 5. The effect of anisotropic pitch angle diffusion on the spectral index produced 
by a parallel relativistic shock front. The full line (isotropic) and the short dashed line 
(anisotropic) show the results of the Qj computation J69[ , with, for the latter, pitch 
angle diffusion given by equation (52). The spectral index found by Heavens and Drury 
[ [70[ is shown by the dotted line (isotropic) and long dashed line (anisotropic). The 
pitch angle diffusion coefficient in the latter case is given by equation (53). Note that 
[|39| and j70| assume slightly different compositions for the plasma. 
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3.3.1. Simulation of magnetic field The transport of charged particles in a relativistic 
flow can be simulated in a number of different ways. At perhaps the most basic level, 
the entire electromagnetic field can be specified in the relevant regions of space, and 
the particle orbits integrated directly. The basic idea behind this approach is that 
the waves responsible for the transport of energetic particles are low-frequency, MHD 
modes, with phase velocities much smaller than the velocity of the shock front. In 
this case, it seems a reasonable approximation to neglect any time dependence of the 
electromagnetic field, apart from that associated with the fast-moving discontinuity 
which is the shock front itself. A single realisation of a stochastically generated static 
magnetic field configuration is then sufficient to determine the motion of all particles. 

This method has been employed by two groups to investigate particle acceleration 
at relativistic shocks []73| , [74], with contradictory results. There are several problems 
with the method. The most severe of these is the limited dynamic range available for 
simulation of the stochastic magnetic field. Ballard & Heavens |73| use a 3-dimensional 
turbulent field with vanishing average component, generated by waves of random phase 
and intensity over a wavelength range of a factor of 100. Thus, a particle which during 
acceleration changes its gyro-radius by more than this factor will feel the effects of the 
limited dynamic range of the wave turbulence. Such a small energy range means that 
the characteristic index of the power-law spectrum expected of a realistic system can be 
determined only very approximately. Ostrowski ]74|], on the other hand, uses a larger 
dynamic range, but employs only three monochromatic turbulent waves of stochastic 
amplitude but fixed phase within this range. He presents results for several orientations 
of the average magnetic field with respect to the shock normal. However, it is not clear 
that the acceleration process is independent of the number of wave modes used. 

Inherent in each of these methods is another difficulty: suppose that we succeeded in 
realising a stochastic magnetic field over a sufficiently wide range of length scales. How 
can we be sure that the numerical orbit integration reproduces the particle trajectories 
with sufficient accuracy? At first sight, this does not appear to be a serious problem. 
However, recent studies of anomalous transport regimes in fusion plasmas (e.g., Rax & 
White |75|| ) and in the non- relativistic acceleration problem [76 show that the property 
of reversibility of a particle trajectory in a static magnetic field is both difficult to achieve 
numerically and important for the determination of the characteristic spectral index in 
shock acceleration. The problem can be understood by considering the situation when 
long wavelength fluctuations create a local magnetic configuration in which a particle 
is 'trapped'. When such a structure passes over the shock front, the energy gained 
by a trapped particle depends strongly on how effectively the particle is confined |77[ . 
Physically, one might expect the perturbing effect of non-static field fluctuations to 
play a role. This can be interpreted as an effective 'decorrelation time' over which the 
particle forgets about a particular realisation of the static stochastic field. Numerically, 
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decorrelation is likely to be determined by the finite accuracy with which particle 
paths are computed. Because of the stringent requirements on accuracy, it has proved 
particularly difficult to devise numerical schemes which allow the study of the anomalous 
transport properties which arise on a time-scale shorter than the decorrelation time |75|] . 



3.3.2. Simulation of particle transport An alternative, but perhaps less fundamental, 
approach to numerical simulation of acceleration is to make an assumption not about 
the stochastic field itself, but about the transport properties which it produces in the 
particles which are to be accelerated. This was a crucial step in developing the original 
analytic diffusive shock acceleration models. It results in a kinetic equation governing 
the single particle distribution function which contains a collision operator. 

The basic idea of the Monte Carlo method is to find a way of constructing a 
stochastic trajectory whose distribution obeys the desired transport equation. Then, 
by repeating the procedure a large number of times, the distribution itself can be 
constructed approximately. This approach to solving transport problems has a long 
tradition, particularly in the fields of neutron transport and radiative transport. 
The transport of charged particles moving in a turbulent magnetic field is, however, 
qualitatively different, since the turbulence acts continuously to deflect the particle, 
whereas the scatterings of a neutron or photon are impulsive events which act on the 
trajectory during short time intervals between which are long periods of undisturbed 
flight. The differences are reflected in the mathematical properties of the collision 
operator contained in the transport equation. The pitch-angle diffusion operator 

c ^ ' (54) 

results in a transport equation which is a second order differential equation, for example 
that given in equation (42). This describes the continuous deflection of the particle 
by an infinite succession of infinitesimally small changes in pitch angle. On the other 
hand, the transport of photons and neutrons is described by scattering through a finite 
angle. In the simplest case of elastic scattering, one has, ignoring for the moment the 
gyro-phase, 



C [/(/*)] = ^scatt 



1 



d//p(/i, //)/(//) - f{n) 



(55) 



where /i sca tt is the scattering frequency and p(/i, //)d/z is the probability that in a single 
scattering a particle with pitch angle cosine // is scattered into the range dfi around /i. 

The spectral properties of these two operators are quite different pi] , as is the 
approach to constructing a corresponding stochastic trajectory. In the first case, 
one formally makes use of the Ito calculus (e.g., Gardiner [|7£|) to write a stochastic 
differential equation (SDE). The appropriate stochastic trajectory is then found by 
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numerically integrating the SDE. An intuitive picture of this procedure is the following: 
to proceed from one point on a stochastic trajectory (labelled by, say x^, Hi) to the 
next, one must solve the transport equation under the initial condition f{x,fi,t = 0) = 
S(x — Xi)5(fi — fii). For small changes in x, this can be done approximately by assuming 
that fi changes only slightly. Expanding equation (42) in powers of h ~ Hii one finds 

r ± (l + v±m)?£ + T ± (v ± + ^) d J- = DM-^f + D'^)^f (56) 

where D' = dD/dfi. The substitutions A = (t— ti)/[r±(l+/j li v±)), 5 = (x— Xi)/\T±(v±+ 
Hi)] and r] = [i — Hi + D'A reduce this to the heat conduction equation, and the solution 
is easily seen to be 

5{x -Xi — T±(v± + Hi) A) 



f[x, fi,t = ti + r±(i + iiiV±)A]) 



ttD^A 



exp 



-(H-Hr-D'^A) 2 /(D^A)\(57) 



The next point on the trajectory is found by setting a sufficiently small time step 
A and choosing a new stochastic value of h from the Gaussian distribution of 
equation (57). This method, together with specialised techniques for enhancing the 
statistical significance of the results (such as the 'splitting' technique) was applied to 
the the acceleration problem for particles undergoing synchrotron losses in |79j . The Ito 
approach to more general problems including second order Fermi acceleration has been 
presented by Kriills & Achterberg p]|j . 

Whereas the Fokker-Planck operator stems from the quasi-linear theory of plasma 
turbulence, as discussed in Sect. |372"| , there is currently no theoretical justification for 
a large-angle scattering operator such as equation (55) to describe charged particles 
propagating in a turbulent plasma. Nevertheless, motivated by the qualitative 
appearance of particle trajectories in simulations of strong plasma turbulence ||81|1 , 
which show sudden, large-angle deflections, a mixture of the two operators has been 
investigated analytically [p4jl . The motivation to perform Monte-Carlo simulations 
using the large-angle scattering operator is even stronger, since it is both relatively 
simple to implement and modest in its use of computer time. A stochastic trajectory 
is constructed by first choosing a scattering time distributed exponentially in the usual 
manner: At = — ^scatthi^), where £ is a random number distributed uniformly on the 
interval [0, 1]. During this interval, the particle propagates with unchanged pitch angle 
from position to x i+ i = Xi + (v± + Hi) At. The new value of /i = fa+i is then found 
from the distribution p(fi,fa)- 111 t ne non-relativistic case, this method has proved 
highly popular (for a review, see |52]| ). Despite the lack of a theoretical background, the 
results it yields agree with solutions found using the diffusion approximation. This is 
because the distribution function in the non-relativistic case is close to being isotropic. 
Independent of the mathematical structure of the collision operator, it is then always 
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possible to derive the equation of spatial diffusion, so that a Monte-Carlo simulation 
using equation (54) must yield the same results as one performed using the more 
economical operator equation (55). 

Once the motion becomes relativistic (or, more generally, once the fluid speed 
approaches the particle speed), anisotropics arise, and, with them, differences between 
the pitch-angle diffusion and large-angle scattering cases []64|, |B3[| . Nevertheless, large- 
angle scattering with complete isotropisation at each scattering (i.e., p(/i, //') = 0.5) 
has been used to estimate the acceleration time at relativistic shocks []8J], |83], ^ a 
quantity which is needed in order to estimate the maximum energy cosmic rays which can 
be produced by the relativistic jet of a radio galaxy. It is, of course, possible to choose 
the large angle scattering operator such that not only the pitch angle is affected, but 
also the phase and, more importantly, the position of the guiding centre. Particles can 
then move across field lines as a result of scattering, mimicking the quasi-linear process 
of cross-field diffusion. This effect has been studied in connection with acceleration at 



non- relativistic oblique shocks [ 87 1 . If, in addition, the scattering is chosen to be strongly 
peaked in the forward direction, then pitch-angle diffusion is simulated. In this way, 
a comparison of the effects of pitch angle diffusion and large-angle scattering has been 
made in investigations of relativistic, oblique shocks, both with |60], |56| and without [|59[ 



explicit implementation of the assumption of magnetic moment conservation. 

Monte-Carlo simulations have also been performed for highly relativistic shocks by 
Bednarz & Ostrowski [72| using upstream Lorentz factors up to roughly 240. They find 
the spectral index of accelerated particle converges to the value s = 4.2, independent 
of the orientation of the magnetic field, provided both pitch angle scattering and cross- 
field diffusion are permitted. However, in the ultra-relativistic limit, a particle which 
manages to cross a shock front from the downstream side into the upstream flow is very 
rapidly overtaken again once it is deflected. In fact it, can perform only a small fraction 
(~ l/r_) of a gyration about the magnetic field line, unless the direction of the field 



is exactly along the shock normal j45|. In this case, the upstream transport cannot be 
represented by large-angle scattering. Instead, a combination of motion in a uniform 
field, and diffusion in angle due to fluctuations in the field on length scales much shorter 
than a gyro-radius arises. If the field fluctuates rapidly, one would expect to recover 
the operator equation (54), where now the quantity fi is interpreted not as the cosine 
of the pitch angle, but as the cosine of the angle between the particle velocity and the 
shock normal. Gallant et al [35] have extended the method to the ultra relativistic limit 
(r_ — > oo) and considered both the case of diffusion in angle and scatter-free deflection 
by a uniform field. The corresponding power laws are s = 4.25 and 4.3 respectively, 
in reasonable agreement with Bednarz and Ostrowski, despite differences found in the 
angular distribution of the particles. This result is particularly encouraging for those 
theories of gamma-ray burst sources which use a relativistic blast wave to accelerate 
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the particles, since it is close to the index of the particle spectrum required to produce 



afterglow spectra [gg, ^9 
4. Discussion 

At the outset of this paper we pointed out that recent observations of blazars and 
gamma-ray bursts (GRBs) have each given renewed impetus to the study of particle 
acceleration in relativistic flows. We now address the issue of to what extent the results 
of these studies have been brought into contact with the observations of such objects. 

Shock fronts arise when the relativistic flows encounter the ambient material. 
Variability in the central 'engines' of these sources may also result in internal shock 
waves, for example along the jets of AGNs. The shock jump conditions are then 
contained in the theory of section 0, at least at the MHD level, and, given some 
knowledge or model of the constituent matter and magnetic field, the possible shock 
configurations follow. In principle it is then possible to apply the acceleration theory of 
section [5] to model the energisation of highly relativistic particles. However, although 
reliable computations of flow patterns are now available, the most common approach 
to acceleration is still to assume a power law distribution of particles. In those 
instances where the acceleration is modelled at a deeper level it is usually with a simple 
phenomenological equation for the distribution |K| [JT], [52|, [^J. 

Before discussing the issues relating to jets and GRBs it should be noted that there 
exists in Nature one relativistic shock front which we are particularly well-placed to 



observe ions |25| . It is located in the Crab Nebula, and is responsible for thermalising the 
energy which streams from the central pulsar in the form of a relativistic flow of positrons 
and electrons, possibly contaminated by. It is thought that the shock is stationary and 
the magnetic field oriented perpendicular to the shock normal |26|. Using l| D numerical 



simulations, a determined effort has been made to understand both the plasma physics 
of the shock front and the connection with particle acceleration. However, this source 
seems fundamentally different to gamma-ray blazars and gamma-ray bursts. Not only 
is it much less powerful, but it also does not display the rapid variability characteristic 
of the extra-galactic sources. In view of this, and because of the very different technical 
approach employed, we do not attempt to synthesise this theory into our discussion, but 
refer the interested reader to the original papers ^7], |S|. 



4-1. Relativistic Jets 

The problem of how particles are accelerated to non-thermal energies in relativistic jets 
has been discussed for a number of years |94]. Diffusive acceleration at non-relativistic 



shocks is a possibility which appears to provide a reasonable picture of acceleration 
at several jet hot-spots which emit synchrotron radiation. In particular, the observed 
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softening of the spectrum towards higher frequencies has been compared with analytic 
computations which include energy losses by synchrotron radiation p3, |9^] . However, 



both apparent superluminal motion and the rapid variability of emission from blazar 
jets imply substantial Doppler boosting, which indicates that the non-relativistic theory 
is inadequate. 

The first attempt to model the radio emission of relativistic jets was made by Wilson 
and Scheuer PS |, who performed numerical hydrodynamic simulations in which it was 
assumed that a fixed fraction of the energy thermalised at a shock front is converted 
into relativistic particles and magnetic fields. This is a rather sweeping assumption, 
covering not only particle acceleration, but also a dynamo process to convert part of the 
available energy into magnetic flux. Although such processes may well occur, there is 
currently no theoretical prospect of checking whether or not the assumed values of the 
parameters are realistic. Subsequent investigations have chosen to inject particles at a 
particular point and simply allow them to be compressed adiabatically at shock fronts 
according to Eq. (39) |99|, pL QQ| , |101|| . Since all of these are hydrodynamical schemes, 



which do not follow the evolution of the magnetic field, equipartition was assumed. If 
correct at injection, this assumption continues to hold through homologous compressions 
or expansions, but it is violated at a shock front. Of course, if the field is in some sense 
random, the average effect of a shock may be similar to homologous compression, but this 
is not going to be the case for sources with a high degree of polarisation. Only a rapidly 
acting dynamo (or, in some cases, reconnection) could then restore equipartition and 
thus rescue the simulations. Dropping the assumption of equipartition by assuming the 



fields to be frozen into the plasma (as in ideal MHD), Matthews and Scheuer ||102|, |1 03 



have computed the emission patterns, once again assuming that relativistic electrons are 
compressed adiabatically at shock fronts. Although they enable interesting simulated 
'radio maps' to be constructed, the major shortcoming of these simulations is that they 
do not include enough of the physics of particle acceleration to be able to predict spectra. 

The best developed model of blazar spectra is based on the picture of a shock front 
moving down a jet ||104j| . Both the magnetic field and the particle distribution vary 



with position and the observed radiation is a superposition from different parts of the 



jet [105, 106 1 . These models are generally called 'inhomogeneous'. However, recent 



results concerning the rapid variability of blazar sources at all frequencies [|107|| and 



the detection of their emission at energies up to at least 10 TeV [[108|| have provided 
new restrictions on the possible acceleration mechanisms. The observed simultaneous 
variations in X-rays and TeV gamma-rays indicate that a single population of particles 
in a relatively localised region is responsible. Consequently 'homogeneous' models have 



been widely discussed [109, |110|, 111, 112, p3| and there is rough agreement on the 



parameters of the emission region. Furthermore, simultaneous variability of emission at 
widely differing wavelengths is more easily interpreted in terms of directly accelerated 
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electrons than in hadronic models ||113|| , where the high energy emission results from 
the acceleration of protons. 

In the case of blazars the synchrotron emission in the radio to infra-red range 
typically displays a very hard power-law (a < 0.5). For homogeneous models this 
means a flatter particle spectrum is needed than the canonical s = 4 result of non- 
relativistic theory, since the cooling time of the responsible particles is long, and the 
indices are related by a = (s — 3)/2. For Markarian 421, for example, a value of s = 3.7 
is required to give the correct index of a ~ 0.35 for the synchrotron emission. From 
Fig. ^|, we see that, in the case of parallel shocks, such values are obtained only for mildly 
relativistic speeds in plasmas in which the electron pressure dominates. A perhaps more 
plausible alternative is to appeal to the fact that oblique relativistic shocks, through a 
combination of shock-drift and first-order Fermi acceleration, give much harder spectra 



ranging up to s = 3 22 1. This is due to the increased importance of reflections back into 



the upstream region from the shock front with increasing obliquity and shock speed, as 



discussed in section 3.2 



4-2. Fireball models of GRBs 



The most popular, generic model for GRB sources is that of the fireball [[114] , |115| , |116 
which involves the sudden release of energy in the form of an optically thick e^ plasma. 
While initially confined to a small volume, the plasma, and along with it the radiation, 
expands out to the point where the system becomes optically thin, releasing the 
electromagnetic radiation. The presence of a baryonic component can, however, result 
in the conversion of the initial internal energy into the kinetic energy of a thin shell 
of matter which expands relativistically, if such baryonic 'contamination' is small, into 
the surrounding medium ||117|| . When this shell interacts with the ambient medium a 
relativistic shock front forms which ultimately dissipates the kinetic energy of the ejecta 
into post-shock thermal energy and the production of energetic particles. In qualitative 
terms at least this picture is analogous to that of an expanding supernova remnant 
where the theory of diffusive shock acceleration can be applied. However, although the 



hydrodynamics of an expanding fireball have been well studied [14, 11S, 119 , 120 , 121 ], 
the theory of particle acceleration at relativistic shocks, with the exception of some 
Monte-Carlo studies, has not been applied to these sources. This is also true of the 
case where a fluctuating Lorentz factor for the ejecta leads to the formation of internal 
shocks [|122j , |123|| . Likewise, arguments favouring a common origin for ultra-high-energy 
cosmic rays and GRBs ||124| , |125|| have used simple estimates of the acceleration processes 
which may take place. If the problem of acceleration is essentially treated as one 
of injection where a value of s = 4 is assumed and the total energy in energetic 
particles is parameterised as some fraction of the total fireball energy, estimates of 



34 



the synchrotron and synchrotron self-Compton fluxes can be made ||126| , |127|| and 
compared with observations. However, there are many important physical questions 
still unanswered. How does the environment and its equation of state influence the 
spectrum? Can a description of acceleration at multiple, internal shock waves account 
for the variability in GRBs? 

It is not hard to see why the acceleration part of the problem has been largely 
simplified in the literature for GRBs. Until recently, published results were restricted 
to Lorentz factors of less than about 10, whereas GRBs require values ten times larger. 
First results have now appeared for the first-order Fermi process in the ultra-relativistic 
case ]72|. |45|] and report a characteristic value of the power-law index encouragingly close 



to that which the observations seem to indicate gj| ||, |128| . As in the non-relativistic 
case, the compression ratio at an ultra-relativistic shock is independent of the shock 
speed. If, as the simulations suggest, the power-law index of accelerated particles turns 
out also to be independent of the shock speed, modelling of individual events will become 
simpler, and the theory of particle acceleration may play an important role in advancing 
our understanding these sources. 
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Appendix A. Hydrodynamic jump conditions 

Consider a fluid made up of a number of constituents (electrons, protons and atoms of 
H and He, for example), labelled by the suffix i, and assign to each a temperature Tj. 
For an ideal gas of monatomic, non-degenerate constituents, the contribution Wi of each 
of them to the total enthalpy density is 

Wi = miUiG^i/Ti) (Al) 

where rrii is the rest-mass of the i'th constituent and rii its proper number density. Note 
that we set the speed of light and Boltzmann's constant equal to unity k-Q = c = 1. The 
function G(z) is given in terms of modified Bessel functions of the first kind (McDonald 
functions): 

g « - W, (A2) 



35 



129||), and has the asymptotic expansion for low temperatures 



5T 

GimlT) -»• 1 + — as T -> (A3) 
2m 

and for high temperatures 

4T 

Gim/T) -> — + 1 as T -> oo (A4) 
m 

Each constituent satisfies the ideal gas law 

^ = riiTi (A5) 
and the overall pressure, enthalpy density and energy density are simply 

P = J^Pi 

W = 

e = w — p (A6) 
We further define the rest-mass density p 

P= Yl m i n i ( A7 ) 

and the abundance (by mass) of the constituent v. 

7] i = m i n i /p (A8) 
so that J2Vi = 1- 

The equation of state for this dissipation free (adiabatic) case is then 

w = p[j2ViG(m i /T i }] (A9) 

together with a relationship between the Tj. (For example, in thermodynamic 
equilibrium: Tj = T, or, for cold electrons/hot ions: T e = 0, T H = T He = T.) 

The speed of sound is simply given by = (dp/de) s where s is the specific entropy 
defined by the first law of thermodynamics: Tds = d(e/p) +pd(l/p). Using Eq. (3), 
and assuming 7 to be constant, we obtain 

(7 - l)p 



2 IP 

"• = 7 



(7 - l)p + 7^^ 



(A10) 



which yields the familiar non-relativistic result for p p and gives v s — > 1/ for a 
relativistic gas. 

The jump conditions Eqs. (7) and (8) can be simplified by the introduction of the 
following notation: let the quantities <fi± be defined by 

r± = cosh0 ± . (All) 

Then one has 



36 



v± = tanh0 ± (A12) 

The possibility of particle creation at the shock front may be taken into account by 
reserving the notation p + for the density of conserved particles, and introducing the 
parameter 77 to relate the total downstream rest-mass density ptotaR t° P+ : 

V = P+/Ptotal+ (A13) 

If particle number is conserved, rj — 1, for small rj, a large number of particles of non- 
zero rest-mass are created at the shock. In terms of the energy per unit rest mass 
e_ = e_/p_ and e + = e + /p to tai+ and the similarly defined pressure per unit rest mass 
P±, one obtains 

p + sinh</> + = p_sinh0_ (A14) 

(e+ + P + )sinh 2 0+ + P+ = V~ \(e- + P-)sinh 2 0- + Pi (A15) 

P+ L J 

(e+ + P + )sinh</> + cosh0_|_ = 77— sinh0_ cosh0_(e_ + PL). (A16) 

P+ 

Elimination of p + / p_ leads, after a little manipulation, to 

u; + cosh0 + = r]W^cosh.(f>- (A17) 



and 

(w + cosh 2 (f) + — e + ) = ?7e_sinh0 + sinh0_ 1 H — — I , (A18) 



P_coth 2 



where w± = e± + P±. The last term on the right-hand side of equation (A18) can be 
neglected if one restricts one's consideration to strong shocks, i.e. to shocks in which the 
upstream pressure is negligible. This approximation may be extended by assuming the 
upstream medium is cold, in which case the gas possesses only that energy attributable 
to its rest mass: e_ = 1. Then, eliminating cosh0 + , one has 

coshv = y+-^ 2 ) (A19) 

7/ 2 (e+ — T] 2 — Pf) 

which is a generalisation to the case i| ^ 1 of a relation found by Peacock |57|| . 
Alternatively, one may obtain from Eqs. (A17) and (A18) the relation 

e + = r/cosh0_cosh0 + — ^sinh0_sinh0 + 

= 77cosh(0_ - 0+) (A20) 

The right-hand side of Eq. (A20) is just 7] times the Lorentz factor r re i associated with 
the relative velocity of the upstream fluid with respect to the downstream fluid. Thus, 
in the absence of particle creation (77 = 1), the average energy per particle as seen from 
the frame of the downstream fluid is constant across a strong shock (in which e_ = 1). 
ODD. 
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For a given upstream velocity t>_, equation (A19) must be solved numerically to 
give the downstream parameters. The procedure is as follows: v _ determines cosh0_ 
through Eq. (All). The right-hand side of (A19) is a function of the single parameter 
e + (the average Lorentz factor), since the pressure P + is given through the equation of 
state Eq. (A9). A root finding algorithm can therefore be used to find e+. By means of 
equation (A17), which for strong shocks specialises to 

u; + cosh0 + = ?7cosh0_, (A21) 

one finds the downstream velocity v + . 

The ultra-relativistic case is best considered without imposing the restriction to 
strong shocks. The equation of state p = e/3 can then be applied in both the upstream 
and downstream regions. Returning to equations (A17) and (A18), one finds, after a 
straightforward calculation, 

v-v + = i (A22) 

A useful analytic expression for strong shocks can be found if the equation of state 
(3) is used [14|]. The procedure is to write the specific enthalpy in terms of r rel using 
Eq. (A20): 

w + = 7(r?r rcl - 1) + 1 (A23) 

Remembering that W- = 1 for a strong shock, the cosh0 + term of Eq. (A17) can be 
rewritten as a function of <p + — 0_ and </>_ to give an expression for tanh0_. Writing 
this in terms of r re i = cosh(0 + — 0_) one finds 

- 7(2-7)(r/r rel -l) 2 + 2( V T rel - 1) + 1 - r/ 2 1 ; 

Thus, given r] and the relative speeds of the upstream and downstream fluids, the shock 
speed and downstream pressure can be written down. 



Appendix B. MHD jump conditions 

In solving the MHD shock jump conditions it is useful to define the components of b± 
and u± along the shock normal, d± = b±l^ and a± = u±i M . These Lorentz scalars take 
the values d± = T±(B±/ ^AnxpS) cos$± and a± = T±v±. With n± = p±/p~ the jump 
conditions become 



[na\ 

[v a ] 

[W a 







[du a - ab a ] = 

n / 1 

anGu a + -l a + \b\ 2 au a + -l a 
z V 2 



-db c 



(B25) 
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where z = p/p. These are a set of non-linear equations in the unknown downstream 
quantities. However, they can be reduced to one transcendental equation for z + , which 
is solved numerically, from which the other downstream quantities can be deduced. The 
first step in obtaining the equation for z + involves the invariance across the shock of 
na, V a V a , W a V a , W a l a - V a V a and X a X a where X a = W a - (a 2 nG + n/ z + \b\ 2 /2)l a . 
This leads to five equations for the five unknowns a+, 



l&R , d + , n + and z + . These are 



na 
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a 2 \b\ 2 
[dG] 
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d 2 



, r . n \b\ 2 

z 2 



2 2 

a n 



G 2 (l + 
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z v 



a 2 \b\ 2 ) + 2a 2 \b\ 2 nG 













(B26) 
(B27) 
(B28) 

(B29) 
(B30) 



For algebraic simplicity we introduce the quantities M = na, Q = (a 2 |6| 2 — d 2 )/M, 
P = dG/ \[M and D + = (Q + P 2 /G 2 + )/2 and the five equations can be written as [plj 



[M] = [Q] = [P] =0 
1 D + 

a + G + H 1 2~ = ci 

Cb-Sf-Z-Sf- 0_|_ 

(1 + a 2 + )G 2 + + — Ud + G + -9\ = C2 
a + \ z+) 



(B31) 
(B32) 

(B33) 



where c\ and C2 are determined by the given upstream state. Equations (B32) and (B33) 
contain just the two unknowns a + and z + since D + and G + are functions of z + and all 
other quantities are given by the upstream parameters. These two jump conditions can 
be manipulated into a linear equation for a + 



(c 2 + G + /z + - G 2 + ) p - (a + D + G + ) ac t 



c\(3G + — a 2 G A 



(B34) 



where 



a = 3D+G 



2Q 



and P = 



AD + G^ 



2Q 



D+G 2 + 



D+c 2 . 



Inserting into either (B32) or (B33) gives a transcendental equation for z + which can 
be solved by a root finding algorithm. 

For a dynamically unimportant magnetic field the flow velocities are directed along 
the shock normal so that Eq.(28) holds. A second relation between upstream and 
downstream magnetic fields can be obtained from (B27) which becomes 



T 2 B 2 



cos 2 <3>_ 



T 2 + B 2 + (t 



cos 2 $_ 
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When combined with Eq.(28) this gives the relationship (Eq.29) between B + , the 
upstream parameters and the shock compression ratio which is obtained from the 
hydrodynamic jump conditions. 



Appendix C. The eigenvalues of the relativistic transport equation 

In the diffusion approximation, it is assumed that the distribution function is everywhere 
close to isotropy. Applying this idea to the region far upstream of a shock front suggests 
that some information about the eigenvalue A| =1 and its eigenfunction Q" =1 (//) for small 
v. may be gained by a perturbative approach. This is indeed the case, and one finds 
(assuming is an even function of /i) 

A i f +1 j > 1 - A 4 ' 2 



Q°M = l--±f dfi'—^-. (C35) 

l J-l Lfu'u' 



and 



A^ = 81 



+1 ft ,^2 



dfi — — - 

i D,, 



(C36) 

'mm 

This can be expressed in terms of the spatial diffusion coefficient k, which relates the 
gradient of the particle density to the flux: 

k = c 2 /A\ (C37) 

(where we have explicitly reintroduced the particle speed c) in agreement with the 
expression given, for example, by Skilling | ]131| |). 

In the ultrarelativistic case, particles which cross the shock from downstream to 
upstream have their pitch angles concentrated in a very narrow cone — 1 < // < t>_, as 
seen in the upstream rest frame. This indicates that a change of independent variable 
is appropriate. Choosing y = (1 — v), and approximating D m by the first term 

in its Taylor expansion about n — — 1: w d(l + fi) transforms Eq. (44) into 

yQ" + Q' + e 2 A(y - l)Q/d = (C38) 

(||62|) where a prime indicates differentiation with respect to y and we have simplified 
the notation for Q and A. We have also introduced the small parameter e = 1 — v ~ 
l/(2r 2 ). For large y, the solutions to this equation have the asymptotic dependence 
Q ~ exp(±eyy / A/(i), whereas the regular solution close to the singular point y = has 
the dependence Q oc 1 — (Ae 2 /3d)y + 0(y 2 ). Therefore, bounded solutions are of the 
form 

N 



Q(y) = exp(-ey y/A/d) E a nV\ (C39) 



n=0 
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i.e., a polynomial with an exponential factor. Inserting this into Eq. (C38) and requiring 
that a^v+i vanish yields the positive eigenvalues: 

« d(2i - lf/e 2 . (C40) 



a n +i = jz , T\2 a ™- K^ 1 ) 



In addition, the following recursion relation is obtained: 

2(2N+l)(N-n) 
(n + 1) 

The eigenfunctions Q\{^) are found from these equations by setting N = i — 1 and 
y = (l + n)/(l-v). 
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